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Abstract 


This thesis presents a numerical model that has been developed for determination 
of thermal behavior of a solar collector. The model takes into account the multi- 
dimension and transient aspects that characterise the phenomena of heat transfer 
in a collector. The analysis of different elements of a flat plate solar collector, viz., 
cover, absorber plate, insulation and fluid has been coupled in a global algorithm. 
The algorithm is presented, herein, in the form of a flow chart. The analysis allows 
the prediction of instantaneous temperature distribution in various components of 
a collector while fluid is continuously flowing through it, as well as the performance 
and the influence of various parameters such as mass flow rates, inlet fluid tempera- 
ture and plate tilt. However, the same analysis can be used to predict the influence 
of other parameters such as collector dimensions, material properties etc., also on 
the performance. 
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Nomenclature 


A area (m 2 ) 

abs absorbed 

C specific heat (J.kg' 1 .K~ ] ) at constant pressure 
D{ inner diameter of duct (m) 

D 0 outer diameter of duct (m) 

h heat tranfer coefficient .K -1 ) 

i,j,k indicator of i,j and k th spatial position in x,y and z directions respectively 
k thermal conductivity (W.m~ l .K -1 ) 

I 0 solar radiation flux on the collector per unit area (W.m~ 2 ) 

m total fluid mass flow flux (kg.s -1 ) 

n number of nodal points 

N number of ducts 

q rate of heat transfer per unit surface (W.rn -2 ) 

Q rate of heat transfer (W) 

R perimeter 

S amount of heat flux absorbed by the absorber plate (W.??i -2 ) 

t time (s) 

x,y,z spatial coordinates 
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Greek Symbols 

a absorptivity of absober plate 
ft collector tilt (radian) 

-> surface azimuthal angle (radian) 

(5 thickness (m) 
f einissivity 

0 latitude (radian) 

p reflectivity, density (kg.m -3 ) 

9 angle of incidence of radiation on a tilted surface 
0 2 refracation angle (radian) 

9 Z angle of incidence of radiation on a horizontal surface 
r transmissivity of glass cover 

a Stefan-Boltzman constant (= 5.67*10~ 8 W .m~ 2 .K~ 4 ) 
7] efficiency 

Subscript 

a ambient 
b beam 
c cover 
d diffuse 
f fluid 

1 insulation 
g air gap 

in inlet conditions 
m mean 

o output conditions 
p absorber plate 
w wind 
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Chapter 1 
Introduction 


1.1 Importance of solar energy 

While conventional sources of energy can not be neglected, it is imperative to develop 
programmes in the area of renewable sources as they not only help conserve scarce 
conventional sources of energy, but also contribute immensely to the improvement 
of environment, employment generation, upgradation of health and hygiene, etc. 
The potential for energy generation from these sources -sun, wind , biomass, etc - 
is immense. 

Taking a global energy view, it is felt that solar energy is one of the most promis- 
ing alternate sources because of its diffuse nature and world wide availability. Solar 
energy has a great future in the tropical country like India, where many places re- 
ceive a peak solar flux of the order of 1 kW/m ?. The need for solar energy has 
become especially relevant in the wake of the energy crisis in the rural sector. In 
India, there are still 1,00,000 villages where darkness is not dispelled by electricity. 
These villages ,due to their remote location, are not likely to have access to electric- 
ity supply from the national grid for a very long time partly because of high initial 
cost of transmission. In the rural areas, domestic energy consumption, especially for 
cooking, constitute a significant part of energy demand. The other uses of energy 
are for lighting, water heating, irrigation, post harvesting operations and village 
industries like brick kilns, flour mills, rice mills and agro processing units. Some of 
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ilit' viable uses of solar energies are as follows : 


1.1.1 Solar water heating 

Many flat plate and concentrating collectors have been designed for their usefulness 
in the water heating. They consists of a box with glass covers, covering a metalic 
absorber plate carrying water in the ducts. The water is heated by solar energy 
absorbed by the blackened absorber plate. 

1.1.2 Solar cooking 

A solar cooker can be easily designed and fabricated. The simplest type is a box 
with glass cover facing the sun. Mirrors can be provided to give some degree of 
concentration winch would bring the food, to be cooked, up to cooking temperatures 
quickly. 

1.1.3 Space heating 

Space heating is generally employed in cold countries where keeping the house warm 
is a big problem. This is normally achieved by using solar air-heaters (active solar 
heating) or by designing the house itself to trap the solar energy (passive solar- 
heating). 

1.1.4 Conversion to electricity 

Electricity can be generated directly from sunlight with the help of solar photovoltaic 
(SPY) technology. The basis is the photoelectric effect. The atoms of certain ma- 
terials such as selenium, possess electrons that are easily knocked out of place by 
light energy. Another metal that is in contact with the electron-emitting metal 
collects the electrons and passes them along into wires in a steady stream, while 
other electrons from the wires flow in to replace them. Thus, an electric current is 
established. 
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1.1.5 Solar distillation 

One of the major problems in many parts of the world is the non-availability of 
fresh water. Solar stills can convert saline or brackish water into fresh potable water 
easily, by the process of evaporation and condensation. 

Some other uses of solar energy worth mentioning are solar refrigeration and 
air conditioning, solar furnaces, solar pumps, heating sewage digestors in sewage 
treatment, solar pond, generation of mechanical power, drying grains or vegetables 
etc. 

1.2 Solar collectors and their main components 

In the last few decades, new solar collectors of different shapes have come up. The 
flat plate collectors have evoked considerable interest because they are the simplest 
among the solar energy collecting devices and have a wide range of potential for 
heating and cooling applications. They usually consists of the following components: 

Absorber plate 

It is a blackened metallic plate which is exposed to sun. It is made of a highly 
conducting material with a view to trap maximum amount of radiation incident 
upon it and transferring it to the fluid flowing through the tubes in contact with 
the absorber plate. 

Cover plate 

It is a glass plate usually placed over the absorber plate to reduce heat losses. It 
allows the short wavelength radiation from the sun to pass through it and impinge on 
the absorber plate. Due to the heating up of the absorber plate long wave radiations 
are produced. These radiations are prevented by the glass cover from passing into 
the atmosphere. Hence it leads to higher collection temperature for fixed heat losses. 
The cover plate may be one or tw r o in number. 
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Fluid 


It is a means to extract the heat absorbed by the absorber plate. This hot fluid 
may directly be used for heating purposes or may transfer its heat to another fluid 
through a heat exchanger. 

Insulation 

The bottom part of the absorber plate is separated from the ambient conditions 
bv means of a layer of insulated material. This insulation prevents the loss of heat 
by the plate to the surroundings. Usually the insulation used are glass wool, wood 
wool, saw dust, thermocoal etc. 


1.3 Background Information 

In most of the current bibliography - Hottel and Woertz (1942), Duffle and Beckman 
(1980) - the thermal behaviour of solar collectors is calculated on the basis of models 
that do not take into account the multidimensional and transient nature of the 

i 

phenomenon. When these aspects are analysed, the governing energy equations 
that characterize the thermal behavior of a solar collector are difficult to solve. 
Dillerent. models with simplified assumptions have appeared in an attempt to assess 
• he transient behavior of a collector(Smifh, 198G). These models, however, do not 
take into account the multidimensional and transient aspects of the problem. 

Brian (19(il) used a finite difference method of high order accuracy for the solu- 
tion of three dimensional transient heat conduction problems. Thomas (1949) used 
an algorithm to solve elliptical problems in linear difference equations over a net- 
work. Consequently some work dealing with numerical simulation of solar collectors 
was carried out. Saito and Utaka(1984) presented numerical solutions for a solar 
water heater. They considered a transient one dimensional model. Chiou (1982) 
presented a two-dimensional steady model to study the nonuniformity in the fluid 
flow distribution. A multidimensional, but quasi steady analysis has been done by 
Agarwal (f.9.%’). 
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The conventional one-dimensional analysis for a flat plate solar collector with 
water as the heating fluid is based on the following assumptions : 

1. Steady-state performance. 

2. The headers cover a small area of. the collector surface and therfore can be 
neglected. 

3. The headers provide uniform flow to the tubes. 

4. No absorption of solar energy by the glass covers in so far as it affects the 
losses from the collector. 

5. One-dimensional heat flow through the covers. 

6. Negligible temperature drop through the glass covers. 

7. One dimensional heat flow through the bottom insulation. 

8. The sky can be considered as a black body. 

9. Temperature gradient around the tubes can be neglected. 

10. Temperature gradient in the direction of flow and between the tubes can be 
treated independently. 

11. Properties are independent of temperature. 

12. Dust and dirt on the collector are negligible. 

13. Shading of the collector plate is negligible. 


1.4 Scope of the Present Study 

The main objective of this study is to present a numerical simulation taking into 
account multidimensional and trasient modelling of the different elements of a con- 
ventional flat plate solar collector , viz., cover, duct, plate etc and their coupling 
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in a global algorithm. The Finite difference technique has been employed for the 
solution of the governing partial differential equations. 

The present work takes into account the multi-dimensional and trasient aspects 
of a flat plate solar collector. Solar energy absorbed by the glass cover is also taken 
into account. .This analysis, therefore, leads to the relaxation of five assumptions, 
viz., l 5t , 4 th , 5 th , 7 th and the 10 t/l . Thus; it presents a more realistic physical model. 

A computer program has been developed to predict the performance of a flat 
plate solar collector based on this model. The present work predicts the perfor- 
mance of a conventional flat plate solar collector for a given day, location, material 
properties and geometry. It also predicts the performance of the collector with 
variation in the following parameters : 

• The time of the day 

• Mass flow rate 

• Plate tilt angle 

• Fluid inlet temperature 

While the effect of variation of one parameter, on the performance, is studied, 
other parameters are kept constant. 

The configuration of flat plate solar collector (with dimensions) used in present 
analysis is shown in Figure 1.1. The collector has nine tubes and they are in line 
with the absorber plate. 
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Chapter 2 


Problem Formulation 

2.1 Solar Radiation 

Total amount of hourly global radiation reaching the tilted surface at any instant is 
given by (Sukhatme, 1984) 


I 0 = IbRb + IdRd + {R + Id)Rr 


( 2 - 1 ) 


where. 


h : hourly beam radiation 

It! : hourly diffuse radiation 

Ri, : tilt factor for beam radiation 

cos 0 

cos B z 

Rd ■ tilt factor for diffuse radiation 

1+COS P 

~ 2 

Rr : tilt factor for reflected radiation 

= 

Now, 


lb — I bn COS 0 z 

According to ASHRAE model, for a clear (cloudless) day 
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hm 

h 


Ac M-^£r c ) 

CI bn 


where, 

l bn : beam radiation in the direction of the rays 

6 : : angle of incidence on a horizontal surface, i.e. zenith angle 

A,B,C are constants whose values have been determined by Threlkald and Jor- 
dan (1958). 

Total amount, of heat flux absorbed by the absorber plate is given by 


5 = I b R b (ra) b + (( I d R d ) + (l b + I d ))Rr(Ta) d . (2. 2) 

where, 

(rG')i, : transmissivity-absorptivity product for constant radiation 

( ra)d '■ transmissivity- absorptivity product for diffuse radiation 

2.2 Governing Equations for the multidimensional 
and Transient Analysis 

This section provides a more appropriate and accurate multidimensional-transient 
approach for the analysis of a liquid flat plate solar collector. The mathematical 
description of the multidimensional-transient phenomena in the different collector 
zones is presented here. Properties of the fluid and the materials have been con- 
sidered to be independent of temperature and invariant with time. Depending on 
the relative values of the different dimensions, thermal energy transfer through the 
solid zones have been considered to be two or three dimensional e.g. the heat ex- 
change through the insulation has been considered to be three dimensional while 
that through the plate has been considered to be two dimensional. Temperature 
variation of the fluid flowing inside the ducts is taken as one dimensional in the flow 
direction. Finally, the free convection in the air gap zone, i.e., between the absorber 
plate and the flat plate collector cover is taken into account by means of an empirical 
expression suggested by Buchberg et al. (1976). 
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2.2.1 Cover 


The cover thickness is so small that it is reasonable to consider a uniform tem- 
perature through it. Hence, we get a two dimensional distribution of temperature 
,T c (x,y), in the cover. The governing equation can be obtained from an energy 
balance in a small volume of thickness, S c , of the cover as, 


4 dx‘ dy 2 


ci , Qp-c Qa-c , ^ 9T C 

_l — 1 — 1_ q abs — p c o c - 


dt 


(2. 3) 


where, 

p r : density of the cover 

C c : specific heat of the cover 

k c : thermal conductivity of the cover 

g p _ c : heat flux from absorber plate to the cover 

= average convective heat flux (from absorber plate to coyer) + 

average radiative heat flux (from absorber plate to cover) 

— h (T — T \ 4- a (Tprn -Tcm 4 ) 

— ' b p-c\ J - pm -*• cm) “I /_[_ 4. J__ii 

V <p e c } 

: heat flux from' ambient to the cover 

= convective heat flux (from ambient to cover) + 
radiative heat flux (from ambient to cover) 

= h'w(T a ~ Tc) + etciTsky 4 — T r 4 ) 

q a b s : heat absorbed by the cover per unit volume 
T 

— J ° S c 


The symbols T a , T pm , T cm , e p , e c , K stand for ambient temperature, mean plate 
temperature, mean cover temperature, emissivity of plate, emissivity of cover and ex- 
tinction coefficient respectively. The symbol A p _ c stands for convective heat transfer 
coefficient between plate and cover and is obtained by using the empirical expression 
suggested by Buchberg et al. (1976). An empirical expression for convective heat 
transfer coefficient between plate and ambient, h w , has been proposed by McAdams 
(1954). The heat absorbed by cover, q abs , is calculated using Bouger’s law for a 
partially transparent medium (Duffie and Beckman, 1980). 
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It may bo noted that, the thermal radiation term from ambient to cover is calcu- 
lated considering the ambient air as a black body which emits at sky temperature . 

The sky temperature is calculated from empirical relation given below (Sukhatmc, 1984); 


T Hky - Tam b - 6 

Adiabatic condition is considered on the boundary, i.e. 


dTc 

dn 


= 0 


(2. 4) 


2.2.2 Plate 

As in the case of the cover, here also we consider a two dimensional variation in the 
temperature of the absorber plate i.e. T p (x,y). The governing energy equation can 
be obtained, as before, by a heat balance on an element of thickness, <5 P , of the plate 
as, 






dx 2 dy 




, S _ s~, dip 

S T + 6, - Pl, ° r dt 


where, 

p v : density of the absorber plate 

C v : specific heat of the absorber plate 

: thermal conductivity of the absorber plate 

S : solar flux absorbed by the absorber plate 

tjf.p : heat Hux from fluid to absorber plate 

= h,(T,-T r ) 

Qi-p : ' heat flux from insulation to absorber plate 

— K’l Q z 1 2 = 6{ 

f/ r _ 7 , : heat 11 ux from cover to absorber plate 

/, t r p r r \ i 7^'*'" 1 

-Lpm) ' t J_ « _L Z n 


(2. 5) 
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Heat transfer coefficient between plate and fluid, hf, is calculated using the 
assumption that heat transfer from plate to fluid takes place under ‘constant heat 
flux 1 condition (Wong, 1977). 

Adiabatic condition is considered on the boundary, i.e. 

BT 

af = 0 < 2 - 6 > 


2.2.3 Insulation 

The temperature distribution Ti(x,y,z) in the insulation is considered to be three 
dimensional. The governing energy equation is written as : 


, f d 2 Ti , d 2 T { d 2 Ti s ^ BTi 
k,{ Ox 7 by' dz 2 * fl,c ' at 


(2. 7) 


■where, 

Pi 

Q 
k, 


density of the insulation 
specific heat of the insulation 

thermal conductivity of the insulation 
There is no heat generation or absorption inside the insulation. The boundary 

conditions are: 

• Faces exposed to the atmosphere 


-ki(~) = h w (Ti-T a ) (2. 8) 

• Face in contact with the absorber plate 

Ti = T p (2. 9) 

2.2.4 Fluid 

The fluid flow in the duct is considered to be one dimensional i.e. average values of 
temperature, T/(y), are used in each cross sectional area of the fluid stream. The 
governing energy equation of the fluid flow is: 


12 



( 2 . 10 ) 


where, 

Pi 

Cf 

fn 

N 

<b-J 

R 

A/ 


Qp-f-R 


™ r dTj_ 
N Cf Oy 


= pjAfCf 


dTf 

Ot 


density of the fluid 
specific heat of the fluid 
mass flow rate of the fluid 
number of ducts 

convective heat flux from the absorber plate to the fluid 
h;(T r - T,) 

nDi 


(!)A 2 


The inlet boundary condition for the fluid is 
at y = 0 Ij = Tfi 


2.3 Efficiency. of a Flat Plate Solar Collector 

I'lflicicncy is an important parameter to estimate the performance of a collector. The 
instantaneous cllieiency of a Hat plate solar collector is defined as the rate of useful 
heat gain and the flux incident on the collector. The mean efficiency of collector is 
defined as the ratio of total useful heat gain and the total incident radiation on the 
collector over the period of its operation. 



Chapter 3 


Numerical Modelling 


In this chapter, the methods and the discretized equations used in different collector 
zones are presented. These zones are the cover, insulation, fluid and the absorber 
plate. The discretized elements of cover, plate, fluid and insulation are shown in 
Figure 3.1. The grids are so chosen that corresponding ( i,j) th coordiantes are co- 
incident. The global algorithm for the coupling solution is indicated by a flowchart 
(Fig3.2). Implicit finite difference method is used for the solution of the various 
equations. In the method, the entire region is divided into small grids of finite size. 
Temperatures are calculated on all the nodal points at each time step as we move 
forward in time domain. In the nomenclature used, i,j and k represent parameters 
in x,y and z directions respectively and n is the parameter in time domain. Partial 
derivatives are expressed in the following forms: 


dT 

dx 
d 2 T 

dx 2 

dT 

dt 


Tj+i — Tj-\ 

2Ax 

Tj+1 ~ 2T { + Tj-1 

Ax 2 

rjpn + 1 y-'ri 

At 


.14 



T(i,j- 1) 


(a) Discrete elements in the cover 





Tfl.j- !) 


(b) Discrete clement in the plate 



a y 


(c) Discrete element in the fluid 


;T{i.j,k+ 1) 



• T|(i ,j, k - 1) 

(d) Discrete element in the insulation 

Figure 3.1 : Discrete elements of different components 




3.1 Cover 


A set of discretized equations is obtained by numerical approximation of the gov- 
erning equation(2.3) and its boundary condition(2.4) into each one of the control 
volumes of the grid. The evaluation of the heat conducted is made using Alternative 
Direction Implicit (ADI) technique. The convection heat exchanged with outside is 
taken implicitly and remaining heat exchanges are taken explicitly. 

The Boundary Conditions in discretized form can be expressed as 

at 

i = - 1 T c {i-l,j) = T e (i + lJ) 

j = 1 T c (i,j — 1) = T c (i,j + 1) 

* =. n x T c (i + l,j) = T c (i — l,j) 

j = n y T c (i,j + 1) = T c (i,j — 1) 

Two equations are established for each point of the grid, each one implicit in one 
direction. 

(1) First equation in which the conduction is evaluated implicitly in the x direc- 
tion: 


P T c {i - l,j , n + 1/2) + qT c (iJ,n + 1/2) + rT c (i + l,j,n+ 1/2) = s (3. 1) 


where. 


p = k r Ay- 


q = ~[2k c Ay-^~ + 2p c C c AxAy^~ + h w AxAy] 
Ax At 

r = k r Ay-^— 


s = -k c Ax—[T c (i,j - l,n) - 2T c {i,j,n) + T c (i,j + l,n)] — AxAy[h p - c 

Ay 

(Tp'm - Tcm) + j^^~ ] ~ AxA y[h w T amb + ae c (T? ky - T?{i,j,n))) 

''Cp ' € c ' 

-q ab3 AxAy5 c - 2p c C c AxAy6 c — 
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(2) second equation in which the conduction is evaluated implicitly in the y 
direction: 


pT e {i, j — l,n + 1) + qT c (i, j,n+ 1) +rT c (t, j + l,n+ 1) = s (3.2) 


where, 


p = k c Ax 


Ay 


8 8 

q = -[2k c Ax-^ + 2p c C c AxAy-^ + h w AxAy ] 


k r Ax 


Ay 


s = -k c Ay-~-[T c (i-l,j,n + l/2) -2T c {i,j,n + l/2) + T c (i + l,j,n+l/2)\ 

a(T A - T 4 ) 

AiAy[Vc(Iim - T„) + r VTl Ti 1 “ AxA vK T ° + 

^p + fc ' 


-7 r 4 (?:, j, n + 1/2))] - q abs AxAy5 c - 2p c C c AxAy6 ( 


T e (z, j,n+ 1/2) 
At 


The coefficient matrices of the resultant equations 3.1 and 3.2 are tridiagonal 
and, therefore, can be solved by a TYidiagonal Matrix Algorithm (Thomas, 1994 ). 


3.2 Plate 

The governing energy equation for absorber plate is (2.5). A discretized equation is 
obtained by numerical approximation of the governing equation, using ADI criteria. 
Heat exchanges on cover side and with fluid and insulation are calculated explicitly. 
The Boundary conditions are 

at 

i = 1 T p (i — l,j) = T p (i + l,j) 

j = 1 T p (i,j — 1) = T p (i,j + 1) 

i = n x T p (i + 1, j) = T p (i — 1, j) 

j = n y T p (i,j + 1) = T p (i,j — 1) 
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Two equations are obtained for each point of the grid, each one implicit in one 
direction. 

(1) First implicit equation in x direction 


pT p {i - 1,7, n + 1/2) + qT p {i,j,n + 1/2) + rT p (i+ l,j,n+ 1/2) = s (3. 3) 


where, 


p = 

q = -[2 k P ^y-~ + 2p p C p AxAy^A] 
r = k^yA 

s = -k p Ax-^[T p (i,j - l,n) - 2T p (i,j,n) + T p {i,j + l,n)] - h f RAy 

a(T 4 -T 4 1 

[Tj(i.j,n) - T p {i, j,n)J + AxAy[h p . c (T pm - T cm ) + ] ” ki 

[Ti(iJ,n z - l,n) - Tz(i, j,n z ,n)] - SAxAy - 2p p C p AxAy5 p Tp ^^ — 


AxAy 

Az 


(2) Second implicit equation in y direction 


pTp{i,j - 1,71 + 1) + qTp(i,; j,n+ 1) + rT p {i,j + l,n + 1) = s (3.4) 


where, 


p = k p Ax- ? 


Ay 

q = -[2fc p Ax^ + 2p p C p AxAy^] 

r = k p Ax~~ 

Ay 

s = — Ar p Ay-^A[Tp(i - 1, j, n+ 1/2) - 2T p (i, j, n + 1/2) + T p (i + l,j,n + 1/2)] - hjR 

Ax 
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A»[ T,(i,j,n) - T„(i,j,n )] + AxAy[ - T™) + ‘ffi* , - k,^t 

(r P + ^ ~ J ) A2 

[Ti(i,i,n z - l,n) - Ti(z, j,n z ,n)] - 5 Ax Ay - 2p p C p AxAy(5p- p -- t ’ + ^ 

AAC 

The coefficient matrices of resultant sets of equations are tridiagonal and, there- 
fore, can be solved by a Tridiagonal Matrix Algorithm (TDMA). 

3.3 Insulation 

The discretization is carried out by means of a regular three dimensional grid. For 
each one of the control volume elements, a discretized equation is obtained as a 
result of numerical approximation of the governing equation(2.7) and its boundary 
condition(2.8 and 2.9). ADI method is used for evaluation of the heat conduction. 

The boundary conditions are 

at 

i = 1 Ti(i-lJ,k) = T i (i + lJ,k)-2^[T i (i,j,k)-T amb ] 

j = 1 Ti{i,j - 1, Ac) = T i (iJ + l,k)-2^[T&J,k)-T mb ] 

k = 1 Ti(i,j,k- 1) = T i (i,j,k + l)-2^[T i (iJ,k)-T amb } 

» = n T TiU + hj,k ) = T i (i-l,j,k)-2*tfz[T i (i,j,k)-T amb ] 

J = n y Ti(i,j + l',k) = T i (iJ-l,k)-2^[T i (iJ,k)-T amb ] 

A' = v z k) = T p (Uj ) 

Three equations are obtained for each point of the grid, each one implicit in one 
direction. 

(l)First implicit equation in x direction (move y' forward in time domain) 


pTi(i - l,j,k,n + 1/3) + qTi(i, j, k, n + 1/3) + rT,(i+ 1 ,j,k,n+ 1/3) = 


(3. 5) 


where, 
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p = 


ki 


AyAz 

Ax 


. [2k ^ + 3piCi ^l ] 


ki 


Ax ri 1 At 
AxAy 

Az 

A T A 7 . 

ki — — - 1 ,fc,n) - 2T i (*,j,A: 1 n) + 3i(i,j + l,fc,n)] 


-*» 


Ay 

AxAy 

Az 


[T t (i, j, & - 1, n) - 27} (i, j, fc, n) + 2}(i, j,k+ 1, n)] 


—ZpiC{AxAyAz 


7}(z,j,fc,n) 

At 


(2)Second implicit equation in y direction ( move by another y in time domain) 


pTi(i,j — 1, fe, n + 2/3) + qTi(i,j, k,n + 2/3) + rTi(i,j + 1, A;, n + 2/3) = 5 

(3- 6) 


P 


where, 

, AxAz 
= k r 


Ay 

tn1 AxAz n ^AxAyAz, 
q = -[2^.-^— + 3piCi r? ] 


r = h 


Ay 

AxAy 

Az 


At 


a = - 1, j, k, n + 1/3) - 2 7}(t, j, fc, n + 1/3) + Ti{i + 1J, k, n + 1/3)] 

Ax 

— k i —~—~^r i {i, j,k-l,n+ 1/3) - 2T i (* > j, k,n + 1/3) + !}(», j, fc + 1, n + 1/3)) 
Az 

. a a 7i(i,j,fc,n+l/3) 

-ZpiCiAxAyAz — 


(3) Third implicit equation in z direction ( move by last ^ in time domain) 
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pTi(iJ,k-l,n + l) + qTi{i,j,k,n + l) + rTi(i,j,k + l,n + l) = 5(3.7) 

where, 

AxA y 

AxAyAz. 


P = h 


A 2 


c = -pfc^i+SftCi 


r = h 


Az 

AxAy 

Az 

AyAz 


At 


» = - 1, j, k, « + 2/3) - 2Ti{i, j, k, n + 2/3) + 7-(i + 1, j, A:, n + 2/3)] 


-A:. 


A.x 

AxA. 

Ay 


- l,A;,n + 2/3) - 2Ti{i, j,k,n + 2/3) + T { (i,j + l,fc,n + 2/3)] 


-3piCiAxAyAz 


T(z,j,A:,n + 2/3) 
At 


The matrices corresponding to above sets of equations are tridiagonal. Coupled 
with boundary conditions, the above set of equations is solved one by one using 
Tridiagonal Matrix Algorithm (TDMA). 

3.4 Fluid 

For each one of the finite control volume that result by subdividing the pipe, the 
discretization expiations are obtained by numerical approximation of energy equation 
(2. 10). The internal energy transportation is calculated implicitly. 

The boundary condition is 


at j = 0 Tf(i,j) = ,2). 


T f (i,j,n + 1) 


[mCfTfj/J - l,n + 1) + fi/C/A/Ay^f^ + hjPAyT p (i,j)] 
[mCf + PjCjA + hfPAy] 


(3- 8) 
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Thus we see that solving the above system of equations the temperature distribu- 
tion in different components of flat plate solar collector viz. cover, plate, insulation, 
fluid can be found out at time t + A 1. if we know the temperature distribution at 
time t.The initial temperature distribution is assumed to be known. 

Other Features of the Software 

The climatological data like wind speed and ambient temperatures must be pro- 
vided in the form of data files. The time interval of the data normally of one hour, 
does not necessarily coincide with the one used in the discretization; this makes it 
necessary to work with interpolation techniques. The day has been assumed to be 
cloudless and ASHRAE model is used to postulate hourly beam radiation in the di- 
rection of the rays. The initial temperature of the different elements of the collector 
is taken as the ambient temperature when the calculation is initialized. The accu- 
racy improves when the spatial and temporal discretization densities are increased. 
The numerical solution tends towards an asymptotic value when finer grids are used. 
For the discretization used, further increment in the density of discretization does 
not cause appreciable improvement in the numerical solution. 
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In In data : geometry . collector position and location , material properties . inlet fluid flow condition 
Initial temperature distribution , discretisation grids , climatological data , itarling time , doling time 
mass flow rate. 


ime = Starting time 


Input: Instantaneous ambient temperature, temperature distribution in different 
elements . 


Calculations of solar Insolation and absorbed solar energy by cover and plate . 


Calculation of heat exchange between cover and ambient : empirical expression 
Calculation of heat exchange between plate and cover: empirical expression 


Calculation of cover temperatures : : ADI method 


Calculation of plate temperatures: ADI method 


Calculation of insulation tcmpcriurcs: ADI method 


Calculation of fluid tern pcraturcs: implicit scheme . 


Storing of partial results: cover temperatures, plate temperatures, fluid 
temperatures at all points: instantaneous efficiency 


l ime - I nne 4 A t 



Storing of final results: cover temperatures, plate temperatures, insulation 
temperatures and fluid temperatures at all points: mean efficiency. 


Figure 1. Flow chart, of simulation program 






















Chapter 4 


Results and Discussion 


4.1 Results 

The performance of a flat plate solar collector for a given day, location, position, 
dimensions and material properties depends on inputs, namely, starting time, clos- 
ing time, mass flow rate of fluid, inlet fluid temperature and initial temperature 
distribution in different components of the collector. 

The computations in the present analysis have been done for the 21 st day of the 
sixth month(June) of an year at Kanpur (India). The configuration of the collector 
is same as shown earlier in Figure 1.1. The plate tilt angle for all computations, 
except the one in which effect of tilt variation is observed, is 24°(0.9$) . The study 
bv Morse and Czarnecki (1958) suggests that for above plate tilt the amount of 
annual insolation per unit area is maximum. The starting time and closing times 
have been taken 7:00 hrs and 17:00 hrs respectively. The initial temperature of 
different elements of flat plate solar collector has been assumed to be equal, to that 
of ambient temperature. The dimensions of collector, material properties and other 
parameter values used in computation are given in ’Appendix A’. 

The variation of ambient temperature and that of incident flux on collector 
have been shown in Figures 4.1 and 4.2 respectively. Figure 4.3 shows the results 
of the present anlysis and then comparison with available experimental data. It 
gives the variation of the rise in the fluid temperature with the time of the day. The 



/\moient temperature icemigrauej 



Figure 4.1 : Ambient tei 





figure 4.2 : Incident sol 






experimental values are from Desai (1985) for l 3t of March. To facilitate comparison 
one set of computations have been done with present multidimensional-transient 
model for similar inputs. The pattern of rise in fluid temperature with time is 
similar in both experimental data and that obtained from multidimensional-transient 
analysis. The peak temerature rise predicted by the present analysis and that from 
the experimental data occur almost at the same time, the noon. At noon, the rise 
in fluid temperature, for mass flow rate of 18 kg/h and inlet temperature of 21 °C, 
is 33" C as obtained from multidimensional-transient analysis compared with 28° C 
of experimental data. The discrepancy between the experimental values and those 
obtained from multidimensional-transient analysis could be attributed to various 
factors a few of which are given below : 

(1) the insolation impinging on the plate may be less than the monthly average 
assumed in the present analysis, 

(2) the wind speed may be higher, 

(3) cover and plate materials may not be of very high quality, 

(4) collector may not be well insulated from the sides. 

Figure 4.4 demonstrates the variation of mean plate temperature along the di- 
rection of the fluid flow. The plate temperature increases along the direction of the 
flow because at start, the temperature of fluid is low which increases as it flows 
through the collector. The variation of mean plate temperature along flow direction 
has been plotted for three flow rates, namely, 20 kg/h, 60 kg/h and 100 kg/h . As 
expected, it. is observed that the plate temperatures are higher for lower flow rates. 
This can be explained on account of higher fluid temperatures for lower flow rates. 

Figure 4.5 displays the rise in fluid temperatures with time of the day for three 
different flow rates. The temperture rise first increases till noon and then starts 
decreasing. It is maximum at about 12 noon, almost at the same time when incident 
flux is maximum. For mass flow rates of 20 kg/h, 60 kg/h and 100 kg/h, the 
maximum rise in fluid temperatures are 28 "C, 10°C and 0°C respectively. 

Figure 4.6 shows the variation of instantaneous efficiency with the time of the 
day for the three different flow rates, namely, 20 kg/h, 60 kg/h and 100 kg/h. The 
efficiency rises rapidly from a very small value, at the starting time of 7:00 hrs, 
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Figure 4.4 : Mean plate temperature along flow direction for different flow 
rates 
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Figure 4.5 : Rise in the fluid temperature with the time of the day for different 
flow rates 



Instantaneous efficiency 



Figure d.(> : Instantaneous efficiency with the. time of the day for different 
flow rates 


'll 




till about. ]():()() hrs after which it remains nearly constant till about 15:00 hrs and 
then starts decreasing. This can be explained in following manner. The collector is 
in transient state during start and much of the incident flux is absorbed by cover, 
absorber plate and insulation in raising their temperatures befor they reach in a 
quasi-steady state. Thus the useful energy obtained is low during this period and 
hence the efficiency. From 10:00 hrs till 15:00 hrs, the collector is in a quasi steady 
state and hence efficiency is nearly constant. After 15:00 hrs the collector efficiency 
drops rapidly. The results for the period , 10:00 hrs to 15:00 hrs , are in good 
agreement with the predictions of Agarwal (1996) in his multidimensional-steady 
analysis. But the results for the periods , 7:00 hrs to 10:00 hrs and 15:00 hrs to 17:00 
hrs, diverge from the predictions of Agrawal (1996). This divergence is expected as 
the transient characteristic of heat transfer are dominant during these periods and 
the steady state assumption made by Agrawal (1996) is no longer valid for these 
periods. The pattern of efficiency variation predicted by multidimensional-transient 
is similar to that predicted by Sukhatme (1984). 

Figures 4.7 and 4.8 show the variation of percentage distribution of incident flux 
with the time of the day for mass flow rates of 20 kg/h and 100 kg/h respectively. 
The useful energy curve as expected follows exactly the pattern of efficiency curve 
(see Figure 4.6). The stored energy is positive and sigificant from 7:00 hrs to 10:00 
lirs, implying that transient characteristics are dominant during this period. In the 
afternoon stored energy becomes negative, implying that the collector components, 
viz., cover, absorber plate and insulation are dissipating heat instead of absorbing it. 
From figure 4.8 it is evident, that mass for mass flow rate of 100 kg/h during quasi- 
steady state (i.e. 10:00 hrs to 15:00 hrs) about 75% energy is utilized in heating the 
fluid, 8% is lost from the top of the cover, 2% is lost from insulation faces and the 
optical energy loss is 15%. 

Figure 4.9 shows the variation of mean efficiency of the collector with mass flow 
rate. It is observed that efficiency increases as we increase the flow rate from 20 
kg/h, but becomes almost constant for flow rates higher than 40 kg/h. For very low 
flow rates the fluid is not able to carry the heat as fast as the heat is transferred 
from the plate, this could be one of the reasons for lower efficiency at low flow rates. 
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Instantaneous energy distribution (%) 



Figure 4.7 : Instantaneous energy distribution with the time of the day 



Instantaneous energy distribution (%) 



Figure 4.8 : Instantaneous energy distribution with the time of the day 
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Mean efficiency 



figure 4.9 : Mean efficiency variation with the flow rate 


This further implies that for the optimum utilisation of the flat plate collector the 
mass flow rate should be greater than 40 kg/h. 

As mentioned earlier, a plate tilt of 0.9<£ gives the maximum insolation per unit 
aiea (Morse & Czarnecki, 19o8). Infact an optimum tilt is (<^> — 10°) for the summers 
and (</> + 10°) for the winters are recommended (Kern & Harris, 1975). The latitude 
of Kanpur is about 2 0°. Figures 4.10 and 4.11 show the variation of rise in fluid 
temperatures with the time of the day for mass flow rates of 20 kg/h and 100 kg/h 
respectively. Each of the two figures is plotted for three different plate tilt angles, 
namely, 1G", 24" and36°. It is observed that rise in fluid temperature is higher for the 
case of plate tilt of 1G°(<£ - 10°) and lower for the plate tilt of 36 °(<f> + 10°) when 
compared with the ease for which plate tilt angle is 24°(0.9<^) . The analysis thus, 
further confirms the claims of Kern and Harris (1975) . 

Idle fluid inlet temperature is an important operational parameter which strongly 
influences the peformance of a flat plate solar collector. The effect has been illus- 
trated by plotting the mean collector efficiency with the inlet fluid temperature in 
Figure 4.12. The inlet fluid temperature varies from 32°C to 90°C. It is observed 
that efficiency drops sharply at an increasing rate (from 0.64 to 0.13) with increasing 
inlet fluid temperature (from 32 °C to 90°C ). These results are in good agreement 
with those of Sukhatme (1984). 

4.2 Conclusion 

From the results of the present, study, following inferences may be drawn : 

• The results of present study are within 15 percent of the available experimental 
data. 


• The collector can be assumed to be in quasi-steady state from 10:00 hrs to 
15:00 hrs ,but transient characteristics dominate the performance from 7:00 
hrs to 10:00 hrs and 15:00 hrs to 17:00 hrs. 

• The rise in fluid temperature is maximum at 12 noon 1ST for Kanpur (India). 




Time (lire) 


Figure 4.10 : Rise in fluid temperature with time of the day for difTernet tilt 
angles 




Finni c 4.11: Rise in fluid temperature with time of the day for diflernet tilt 



Mean efficiency 



Figure 4.12 : Mean efficiency variation with inlet fluid temperature 





The collector efficiency remains nearly constant for the period from 10:00 hrs 
to 15:00 hrs. 




• Efficiency of the collector increases with flow rate upto 40 kg/h beyond which 
it remains nearly constant in the laminar flow region. 

• Efficiency of collector decreases with increase in inlet fluid temperature. 

• In summers a tilt angle of <j> + 10° gives the best performance. 

The results discussed in the previous sections are intended to illustrate the pos- 
sibilities and the versatility of the analysis carried out. There are a host of other 
parameters like material properties, day of the year, dimension of the collector etc. 
the effect of which can be studied and optimized with the software developed in the 
present work. The numerical model presented here which simulates a flat plate solar 
collector has been shown to be a powerful tool in the analysis of the thermal perfor- 
mance of a collector. With this software depending on requirement, it is possible to 
make a comprehensive test because it takes into account the multidimensional and 
transitory aspects that characterize a collector. 

4.3 Suggestions for future work 

• Experimental verification of the theoretical results presented in the thesis here. 

• Incorporation of the variation of physical properties with temperature of the 
materials involved. 

• The effect of shadow, number of covers and dust on top of the collector needs 
to be investigated. 

• The similar multidimensional-transient analysis can be taken up for some of 
the novel designs like honey comb collector, double exposure collector and the 
thermal trap collector. 
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Appendix A 


The values of the parameters used are listed below : 

Dimension of the collector 
Width of plate 
Length of Plate 

Thickness of the glass cover (S c ) 

Thickness of absorber plate (S p ) 

Thickness of insulation (Si) 

Dia of pipe (D) 

Material properties of the collector 

Thermal conductivity of the glass cover (k c ) 

Thermal conductivity of the copper absorber plate ( k p ) 
Thermal conductivity of the glass wool insulation (hi) 
Thermal conductivity of the fluid(water) (kj) 

Specific heat of the glass cover ( C c ) 

Specific heat of the absorber plate (C p ) 

Specific heat of the insulation ( C ',) 

Specific heat of the fluid ( Cf ) 

Density of the glass cover ( p c ) 

Density of the absorber plate ( p p ) 

Density of the insulation (pi) 


1 m 

2 m 

2 mm 
2 mm 
7 cm 
1 cm 

1.4 W.m~ 1 .K~ 1 
372 W.m-KK- 1 
0.03 W.m-KK- 1 
0.61 W.m-'.K- 1 
750 J.kg-KK- 1 
419 J.kg-KK - 1 
836 J.kg-\K~ l 
4187 J.kg-KK- 1 
2225 kg.m~ 3 
8300 kg.rrT 3 
36 kg.m~ 3 
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Density of the fluid (p/) 

Einissivity of the carbon black absorber plate (e p ) 
Absorptivity of the carbon black absorber plate (a) 
Einissivity of the glass cover (e c ) 

Discretization grids 

Cover 

Absorber plate 

Insulation 

Fluid 

A/, 

Miscellaneous 

Reflectivity of the surrounding (p) 

Collector plate tilt (f3) 

Latitude of kanpur ((/>) 


1000 kg.m ~ 3 
0.92 
0.92 
0.88 


11*11 

11*11 

11 * 11*8 

11 

GO s 


0.2 

24° 

26.28° 
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